Revised Interoceptive Accuracy Scale (IAS-R)

Analysis

library(tidyverse)
library(easystats)
library(patchwork)
library(lavaan)
library(ggraph)
library(tidySEM)
library(EGAnet)

Study 1 - Reanalysis of Murphy et al. (2020)

Reanalyze their data to confirm the factor structure of the IAS.

Data Preparation

# https://osf.io/3m5nh
df1 <- haven::read_sav("../data/Murphy2020/Study 1.sav") |>
  mutate_all(as.numeric) |>
  mutate(Gender = as.character(ifelse(Gender == 1, "Male", ifelse(Gender == 2, "Female", "Other"))))

# https://osf.io/3m5nh
df2 <- haven::read_sav("../data/Murphy2020/Study 6 IAS.sav") |>
  mutate_all(as.numeric) |>
  mutate(Gender = as.character(ifelse(Gender == 1, "Male", ifelse(Gender == 2, "Female", "Other"))))

# Campos (2022) - Study 1 (https://osf.io/j6ef3)
df3 <- haven::read_sav("../data/Campos2022/Dataset_Test.sav") |>
  mutate_all(as.numeric) |>
  mutate(Gender = as.character(ifelse(Sex == 1, "Male", ifelse(Sex == 0, "Female", "Other")))) |>
  rename(
    Heart = IAS1, Hungry = IAS2, Breathing = IAS3, Thirsty = IAS4,
    Urinate = IAS5, Defecate = IAS6, Taste = IAS7, Vomit = IAS8, Sneeze = IAS9,
    Cough = IAS10, Temp = IAS11, Sex_arousal = IAS12, Wind = IAS13, Burp = IAS14,
    Muscles = IAS15, Bruise = IAS16, Pain = IAS17, Blood_Sugar = IAS18,
    Affective_touch = IAS19, Tickle = IAS20, Itch = IAS21
  )

# https://osf.io/3eztd
df4 <- haven::read_sav("../data/Lin2023/Study 1 & 3.sav") |>
  select(-sex) |> 
  mutate_all(as.numeric) |>
  mutate(Gender = as.character(ifelse(sex_dummy == 1, "Male", ifelse(sex_dummy == 0, "Female", "Other")))) |>
  rename(
    Age = age,
    Heart = Heart, Hungry = HUNGR, Breathing = BREAT, Thirsty = Thirs,
    Urinate = URINA, Defecate = Defec, Taste = TASTE, Vomit = VOMIT, Sneeze = Sneez,
    Cough = COUGH, Temp = TEMPE, Sex_arousal = SEXAR, Wind = WIND, Burp = Burp,
    Muscles = MUSCL, Bruise = Bruis, Pain = PAIN, Blood_Sugar = BloSu,
    Affective_touch = Touch, Itch = ITCH
  ) # No tickle because same Chinese character

# https://osf.io/3eztd
df5 <- haven::read_sav("../data/Lin2023/Study 2.sav") |>
  mutate(Gender = as.character(ifelse(Sex == "男", "Male", ifelse(Sex == "女", "Female", "Other")))) |>
  rename(
    Heart = Heart, Hungry = HUNGR, Breathing = BREAT, Thirsty = Thirs,
    Urinate = URINA, Defecate = Defec, Taste = TASTE, Vomit = VOMIT, Sneeze = Sneez,
    Cough = COUGH, Temp = TEMPE, Sex_arousal = SEXAR, Wind = WIND, Burp = Burp,
    Muscles = MUSCL, Bruise = Bruis, Pain = PAIN, Blood_Sugar = BloSu,
    Affective_touch = Touch, Itch = ITCH
  ) # No tickle because same Chinese character

df6 <- read.csv("https://raw.githubusercontent.com/DominiqueMakowski/PHQ4R/main/study2/data/data.csv") |>
  rename(
    Heart = IAS_1, Hungry = IAS_2, Breathing = IAS_3, Thirsty = IAS_4,
    Urinate = IAS_5, Defecate = IAS_6, Taste = IAS_7, Vomit = IAS_8, Sneeze = IAS_9,
    Cough = IAS_10, Temp = IAS_11, Sex_arousal = IAS_12, Wind = IAS_13, Burp = IAS_14,
    Muscles = IAS_15, Bruise = IAS_16, Pain = IAS_17, Blood_Sugar = IAS_18,
    Affective_touch = IAS_19, Tickle = IAS_20, Itch = IAS_21
  )

df7 <- read.csv("https://raw.githubusercontent.com/RealityBending/IllusionGameReliability/main/data/preprocessed_questionnaires.csv") |>
  rename(
    Gender = Sex, Heart = Item_IAS_Interoception_1, Hungry = Item_IAS_Interoception_2,
    Breathing = Item_IAS_Interoception_3, Thirsty = Item_IAS_Interoception_4,
    Temp = Item_IAS_Interoception_5, Sex_arousal = Item_IAS_Interoception_6,
    Urinate = Item_IAS_Elimination_1, Defecate = Item_IAS_Elimination_2,
    Vomit = Item_IAS_Elimination_3, Wind = Item_IAS_Expulsion_1,
    Burp = Item_IAS_Expulsion_2, Sneeze = Item_IAS_Expulsion_3,
    Muscles = Item_IAS_Nociception_1, Bruise = Item_IAS_Nociception_2,
    Pain = Item_IAS_Nociception_3, Affective_touch = Item_IAS_Skin_1,
    Tickle = Item_IAS_Skin_2, Itch = Item_IAS_Skin_3
  ) |>
  filter(!is.na(Urinate))

Participants

• Sample 1: Data from Murphy’s (2020) study 1, downloaded from OSF, included 451 participants (Mean age = 25.8, SD = 8.4, range: [18, 69]; Gender: 69.4% women, 29.5% men, 1.11% non-binary).
• Sample 2: Data from Murphy’s (2020) study 6, downloaded from OSF, included 375 participants (Mean age = 35.3, SD = 16.9, range: [18, 91]; Gender: 70.1% women, 28.5% men, 1.33% non-binary).
• Sample 3: Data from Campos’ (2022) study 1, downloaded from OSF, included 515 participants (Mean age = 30.7, SD = 10.5, range: [18, 72]; Sex: 0.0% females, 0.0% males, 100.0% other; Gender: 59.6% women, 40.4% men, 0.00% non-binary).
• Sample 4: Data from Lin’s (2023) study 1 and 3, downloaded from OSF, included 1166 participants (Mean age = 32.5, SD = 8.4, range: [16, 60]; Gender: 57.0% women, 43.0% men, 0.00% non-binary).
• Sample 5: Data from Lin’s (2023) study 2, downloaded from OSF, included 500 participants (Mean age = 37.4, SD = 7.4, range: [20, 60]; Sex: 0.0% females, 0.0% males, 100.0% other; Gender: 56.2% women, 43.8% men, 0.00% non-binary).
• Sample 6: New dataset from Makowski (under preparation), included 251 participants (Mean age = 37.5, SD = 13.2, range: [17, 76]; Gender: 53.0% women, 42.6% men, 2.79% non-binary, 1.59% missing; Education: Bachelor, 39.84%; Doctorate, 4.78%; High School, 31.47%; Master, 19.12%; missing, 0.40%; Other, 4.38%).
• Sample 7: New dataset from Makowski (under preparation), included 485 participants (Mean age = 30.1, SD = 10.1, range: [18, 73]; Gender: 50.3% women, 49.7% men, 0.00% non-binary; Education: Bachelor, 45.15%; Doctorate, 1.86%; High school, 34.43%; Master, 15.88%; Other, 2.47%; Prefer not to say, 0.21%).

Total N = 3743.

Descriptive

Distribution

vars <- c(
  "Heart", "Hungry", "Breathing", "Thirsty", "Urinate", "Defecate", "Taste", "Vomit", "Sneeze", "Cough", "Temp",
  "Sex_arousal", "Wind", "Burp", "Muscles", "Bruise", "Pain", "Blood_Sugar", "Affective_touch", "Tickle", "Itch"
)

dens1 <- estimate_density(select(df1, all_of(vars)), method = "kernSmooth") |>
  mutate(Sample = "Sample 1")
dens2 <- estimate_density(select(df2, all_of(vars)), method = "kernSmooth") |>
  mutate(Sample = "Sample 2")
dens3 <- estimate_density(select(df3, all_of(vars)), method = "kernSmooth") |>
  mutate(Sample = "Sample 3")
dens4 <- estimate_density(select(df4, all_of(setdiff(vars, "Tickle"))), method = "kernSmooth") |>
  mutate(Sample = "Sample 4")
dens5 <- estimate_density(select(df5, all_of(setdiff(vars, "Tickle"))), method = "kernSmooth") |>
  mutate(Sample = "Sample 5")
dens6 <- estimate_density(select(df6, all_of(vars)), method = "kernSmooth") |>
  mutate(Sample = "Sample 6")
dens7 <- estimate_density(select(df7, all_of(setdiff(vars, c("Taste", "Cough", "Blood_Sugar")))), method = "kernSmooth") |>
  mutate(Sample = "Sample 7")

rbind(dens1, dens2, dens3, dens4, dens5, dens6, dens7) |>
  ggplot(aes(x = x, y = y, color = Parameter)) +
  geom_line() +
  theme_minimal() +
  labs(title = "Distribution of Responses", x = "Response", y = "Density", color = "Item") +
  facet_wrap(~Sample, scales = "free", nrow = 4)

data1 <- normalize(select(df1, all_of(dens1$Parameter)), range = c(1, 5))
data2 <- normalize(select(df2, all_of(dens2$Parameter)), range = c(1, 5))
data3 <- normalize(select(df3, all_of(dens3$Parameter)), range = c(1, 5))
data4 <- normalize(select(df4, all_of(dens4$Parameter)), range = c(1, 5))
data5 <- normalize(select(df5, all_of(dens5$Parameter)), range = c(1, 5))
data6 <- select(df6, all_of(dens6$Parameter))
data7 <- select(df7, all_of(dens7$Parameter))

data_all <- rbind(
  data1, data2, data3,
  mutate(data4, Tickle = NA), mutate(data5, Tickle = NA),
  data6, mutate(data7, Taste = NA, Cough = NA, Blood_Sugar = NA)
)

Correlations

An overall positive intercorrelation pattern, with no clear structure emerging.

All data

make_correlation <- function(df) {
  correlation::correlation(df, redundant = TRUE) |>
    correlation::cor_sort() |>
    correlation::cor_lower() |>
    mutate(val = paste0(insight::format_value(r), format_p(p, stars_only = TRUE))) |>
    mutate(Parameter2 = fct_rev(Parameter2)) |>
    ggplot(aes(x = Parameter1, y = Parameter2)) +
    geom_tile(aes(fill = r), color = "white") +
    geom_text(aes(label = val), size = 3) +
    labs(title = "Correlation Matrix") +
    scale_fill_metro_c(limit = c(0, 0.75), guide = guide_colourbar(ticks = FALSE)) +
    theme_minimal() +
    theme(
      legend.title = element_blank(),
      axis.title.x = element_blank(),
      axis.title.y = element_blank(),
      axis.text.x = element_text(angle = 45, hjust = 1)
    )
}

make_correlation(data_all)

Sample 1

make_correlation(data1)

Sample 2

make_correlation(data2)

Sample 3

make_correlation(data3)

Sample 4

make_correlation(data4)

Sample 5

make_correlation(data5)

Sample 6

make_correlation(data6)

Sample 7

make_correlation(data7)

Graph Analysis

Exploratory Graph Analysis (EGA) is a recently developed framework for psychometric assessment, that can be used to estimate the number of dimensions in questionnaire data using network estimation methods and community detection algorithms, and assess the stability of dimensions and items using bootstrapping. See https://r-ega.net/articles/ega.html for details.

Unique Variable Analysis (UVA)

Unique Variable Analysis (Christensen, Garrido, & Golino, 2023) uses the weighted topological overlap measure (Nowick et al., 2009) on an estimated network. Values greater than 0.25 are determined to have considerable local dependence (i.e., redundancy) that should be handled (variables with the highest maximum weighted topological overlap to all other variables (other than the one it is redundant with) should be removed).

All data

uva <- EGAnet::UVA(data = data_all, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

 node_i node_j   wto
 Tickle   Itch 0.306

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

  node_i   node_j   wto
  Hungry  Thirsty 0.222
 Urinate Defecate 0.219
  Sneeze    Cough 0.204

uva$keep_remove

$keep
[1] "Itch"

$remove
[1] "Tickle"

Sample 1

uva <- EGAnet::UVA(data = data1, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

 node_i node_j   wto
 Tickle   Itch 0.286
   Wind   Burp 0.275

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

  node_i   node_j   wto
  Sneeze    Cough 0.244
 Urinate Defecate 0.241

uva$keep_remove

NULL

Sample 2

uva <- EGAnet::UVA(data = data2, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

 node_i node_j   wto
 Tickle   Itch 0.350
 Sneeze  Cough 0.317
   Wind   Burp 0.309

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

  node_i   node_j   wto
 Urinate Defecate 0.278

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

uva$keep_remove

$keep
[1] "Cough" "Burp"  "Itch" 

$remove
[1] "Sneeze" "Wind"   "Tickle"

Sample 3

uva <- EGAnet::UVA(data = data3, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

 node_i node_j   wto
 Tickle   Itch 0.445
 Sneeze  Cough 0.318

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

 node_i node_j   wto
   Wind   Burp 0.256

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

  node_i      node_j   wto
  Bruise Blood_Sugar 0.219
 Urinate    Defecate 0.217
   Heart   Breathing 0.208

uva$keep_remove

$keep
[1] "Sneeze" "Itch"  

$remove
[1] "Cough"  "Tickle"

Sample 4

uva <- EGAnet::UVA(data = data4, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

  node_i    node_j   wto
 Urinate  Defecate 0.237
   Heart Breathing 0.218
  Hungry   Thirsty 0.213

uva$keep_remove

NULL

Sample 5

uva <- EGAnet::UVA(data = data5, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

  node_i    node_j   wto
 Urinate  Defecate 0.277
   Heart Breathing 0.267

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

 node_i node_j   wto
   Wind   Burp 0.206

uva$keep_remove

NULL

Sample 6

uva <- EGAnet::UVA(data = data6, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

uva$keep_remove

NULL

Sample 7

uva <- EGAnet::UVA(data = data7, cut.off = 0.3)
uva

Variable pairs with wTO > 0.30 (large-to-very large redundancy)

----

Variable pairs with wTO > 0.25 (moderate-to-large redundancy)

 node_i node_j   wto
 Tickle   Itch 0.297

----

Variable pairs with wTO > 0.20 (small-to-moderate redundancy)

 node_i    node_j   wto
  Heart Breathing 0.209

uva$keep_remove

NULL

Networks

All data

make_egas <- function(data) {
  egas <- list()
  for (model in c("glasso", "TMFG")) {
    for (algo in c("walktrap", "louvain")) {
      for (type in c("ega", "ega.fit")) { # "hierega"
        if (type == "ega.fit" & algo == "louvain") next # Too slow
        egas[[paste0(model, "_", algo, "_", type)]] <- EGAnet::bootEGA(
          data = data,
          seed = 123,
          model = model,
          algorithm = algo,
          EGA.type = type,
          type = "resampling",
          plot.itemStability = FALSE,
          verbose = FALSE
        )
      }
    }
  }

  p <- EGAnet::compare.EGA.plots(
    egas$glasso_walktrap_ega, egas$glasso_walktrap_ega.fit,
    egas$glasso_louvain_ega, egas$TMFG_louvain_ega,
    egas$TMFG_walktrap_ega, egas$TMFG_walktrap_ega.fit,
    labels = c(
      "glasso_walktrap_ega", "glasso_walktrap_ega.fit",
      "glasso_louvain_ega", "TMFG_louvain_ega",
      "TMFG_walktrap_ega", "TMFG_walktrap_ega.fit"
    ),
    rows = 3,
    plot.all = FALSE
  )$all
  list(egas = egas, p = p)
}


egas_all <- make_egas(data_all)
egas_all$p

Sample 1

egas1 <- make_egas(data1)
egas1$p

Sample 2

egas2 <- make_egas(data2)
egas2$p

Sample 3

egas3 <- make_egas(data3)
egas3$p

Sample 4

egas4 <- make_egas(data4)
egas4$p

Sample 5

egas5 <- make_egas(data5)
egas5$p

Sample 6

egas6 <- make_egas(data6)
egas6$p

Sample 7

egas7 <- make_egas(data7)
egas7$p

Structure Stability

Scores

ega <- egas_all$egas$glasso_louvain_ega$EGA
plot(ega)

ega_scores_1 <- EGAnet::net.scores(data1, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

ega_scores_2 <- EGAnet::net.scores(data2, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

ega_scores_3 <- EGAnet::net.scores(data3, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

ega_scores_4 <- EGAnet::net.scores(data4, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

ega_scores_5 <- EGAnet::net.scores(data5, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

ega_scores_6 <- EGAnet::net.scores(data6, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

ega_scores_7 <- EGAnet::net.scores(data7, ega)$scores$std.scores |>
  as.data.frame() |>
  setNames(c("EGA1", "EGA2", "EGA3", "EGA4"))

Factor Analysis

How Many Factors

All data

n <- parameters::n_factors(data_all, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4, 5), ], row.names = FALSE)

n_Factors	Method	Family
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
1	Velicer’s MAP	Velicers_MAP
4	beta	Multiple_regression
4	Scree (SE)	Scree_SE
5	Optimal coordinates	Scree
5	Parallel analysis	Scree
5	Kaiser criterion	Scree

Sample 1

n <- parameters::n_factors(data1, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Bentler	Bentler
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
1	Velicer’s MAP	Velicers_MAP
4	beta	Multiple_regression
4	Optimal coordinates	Scree
4	Parallel analysis	Scree
4	Kaiser criterion	Scree

Sample 2

n <- parameters::n_factors(data2, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Bentler	Bentler
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
4	beta	Multiple_regression
4	Optimal coordinates	Scree
4	Parallel analysis	Scree
4	Kaiser criterion	Scree

Sample 3

n <- parameters::n_factors(data3, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
4	beta	Multiple_regression
4	Optimal coordinates	Scree
4	Parallel analysis	Scree
4	Kaiser criterion	Scree
4	Scree (SE)	Scree_SE
4	BIC	BIC

Sample 4

n <- parameters::n_factors(data4, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
1	Velicer’s MAP	Velicers_MAP
4	beta	Multiple_regression
4	Scree (SE)	Scree_SE

Sample 5

n <- parameters::n_factors(data5, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
1	Velicer’s MAP	Velicers_MAP
4	beta	Multiple_regression
4	BIC	BIC

Sample 6

n <- parameters::n_factors(data6, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Bentler	Bentler
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
1	Velicer’s MAP	Velicers_MAP
4	beta	Multiple_regression
4	Optimal coordinates	Scree
4	Parallel analysis	Scree
4	Kaiser criterion	Scree
4	Scree (SE)	Scree_SE

Sample 7

n <- parameters::n_factors(data7, n_max = 12)

plot(n)

knitr::kable(as.data.frame(n)[n$n_Factors %in% c(1, 4), ], row.names = FALSE)

n_Factors	Method	Family
1	Acceleration factor	Scree
1	Scree (R2)	Scree_SE
1	VSS complexity 1	VSS
1	Velicer’s MAP	Velicers_MAP
4	beta	Multiple_regression
4	Optimal coordinates	Scree
4	Parallel analysis	Scree
4	Kaiser criterion	Scree
4	BIC	BIC
4	BIC (adjusted)	BIC

EFA

All data

efa5 <- parameters::factor_analysis(data_all, n = 5, rotation = "oblimin", sort = TRUE)

Loading required namespace: GPArotation

plot(efa5)

display(efa5)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR3	MR2	MR1	MR5	MR4	Complexity	Uniqueness
Cough	0.72	-5.89e-04	-0.06	-5.43e-03	0.03	1.02	0.50
Sneeze	0.67	-0.03	0.04	0.06	-0.05	1.04	0.53
Burp	0.55	0.13	0.06	-0.03	-0.02	1.15	0.59
Wind	0.48	0.11	-9.06e-03	-1.90e-03	0.06	1.15	0.68
Vomit	0.35	0.03	0.24	0.14	-0.09	2.32	0.68
Itch	0.05	0.71	-0.10	-0.06	0.12	1.12	0.46
Tickle	0.08	0.68	-0.03	0.06	4.30e-03	1.05	0.45
Bruise	8.13e-04	0.45	0.22	0.06	-0.08	1.56	0.68
Blood_Sugar	0.06	0.36	0.14	0.04	-0.06	1.43	0.79
Affective_touch	-0.05	0.36	0.34	0.18	-0.16	2.91	0.64
Muscles	0.03	0.25	0.21	0.17	0.07	2.95	0.72
Taste	0.16	0.18	0.17	0.13	0.05	3.96	0.75
Urinate	0.10	-0.08	0.55	0.01	0.20	1.39	0.55
Defecate	0.22	-0.04	0.46	-0.01	0.15	1.69	0.60
Pain	-0.01	0.22	0.28	0.11	0.23	3.22	0.65
Heart	1.28e-03	-0.03	-0.01	0.66	-0.01	1.01	0.59
Breathing	0.07	7.82e-03	-0.08	0.56	0.14	1.20	0.61
Sex_arousal	0.10	0.07	0.18	0.25	0.12	2.87	0.73
Thirsty	-1.95e-04	0.05	0.11	6.26e-03	0.62	1.07	0.53
Hungry	0.01	0.04	-5.98e-03	0.12	0.60	1.08	0.55
Temp	0.11	0.12	0.14	0.15	0.27	3.12	0.68

The 5 latent factors (oblimin rotation) accounted for 38.36% of the total variance of the original data (MR3 = 10.15%, MR2 = 9.46%, MR1 = 6.49%, MR5 = 6.20%, MR4 = 6.05%).

Sample 1

efa4 <- parameters::factor_analysis(data2, n = 4, rotation = "oblimin", sort = TRUE)
plot(efa4)

display(efa4)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR2	MR3	MR4	Complexity	Uniqueness
Urinate	0.69	-0.18	0.13	-0.04	1.22	0.54
Defecate	0.68	-0.13	0.04	0.06	1.10	0.55
Breathing	0.56	0.03	0.06	7.96e-03	1.03	0.63
Hungry	0.55	0.15	-7.15e-03	0.02	1.16	0.60
Thirsty	0.52	0.05	0.14	-0.01	1.16	0.62
Sex_arousal	0.50	0.14	-0.04	0.20	1.49	0.56
Temp	0.46	0.23	0.16	-0.01	1.77	0.52
Pain	0.46	0.43	-0.05	-0.03	2.03	0.49
Taste	0.40	0.21	-0.07	0.12	1.77	0.69
Heart	0.40	0.06	-0.04	0.15	1.36	0.76
Muscles	0.35	0.24	0.07	0.09	2.07	0.65
Vomit	0.32	0.06	0.23	0.11	2.20	0.68
Tickle	0.02	0.69	0.01	0.04	1.01	0.47
Itch	-0.02	0.67	0.18	0.03	1.16	0.41
Bruise	0.07	0.44	0.15	0.06	1.32	0.65
Blood_Sugar	-0.08	0.38	0.16	0.09	1.61	0.77
Affective_touch	0.38	0.38	-0.13	0.05	2.26	0.63
Sneeze	0.05	0.05	0.71	0.01	1.02	0.41
Cough	0.04	0.02	0.71	0.04	1.01	0.43
Wind	7.52e-03	-0.04	-5.18e-03	0.94	1.00	0.15
Burp	-0.01	0.27	0.17	0.47	1.91	0.48

The 4 latent factors (oblimin rotation) accounted for 44.25% of the total variance of the original data (MR1 = 17.71%, MR2 = 11.72%, MR3 = 7.65%, MR4 = 7.17%).

Sample 2

efa4 <- parameters::factor_analysis(data2, n = 4, rotation = "oblimin", sort = TRUE)
plot(efa4)

display(efa4)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR2	MR3	MR4	Complexity	Uniqueness
Urinate	0.69	-0.18	0.13	-0.04	1.22	0.54
Defecate	0.68	-0.13	0.04	0.06	1.10	0.55
Breathing	0.56	0.03	0.06	7.96e-03	1.03	0.63
Hungry	0.55	0.15	-7.15e-03	0.02	1.16	0.60
Thirsty	0.52	0.05	0.14	-0.01	1.16	0.62
Sex_arousal	0.50	0.14	-0.04	0.20	1.49	0.56
Temp	0.46	0.23	0.16	-0.01	1.77	0.52
Pain	0.46	0.43	-0.05	-0.03	2.03	0.49
Taste	0.40	0.21	-0.07	0.12	1.77	0.69
Heart	0.40	0.06	-0.04	0.15	1.36	0.76
Muscles	0.35	0.24	0.07	0.09	2.07	0.65
Vomit	0.32	0.06	0.23	0.11	2.20	0.68
Tickle	0.02	0.69	0.01	0.04	1.01	0.47
Itch	-0.02	0.67	0.18	0.03	1.16	0.41
Bruise	0.07	0.44	0.15	0.06	1.32	0.65
Blood_Sugar	-0.08	0.38	0.16	0.09	1.61	0.77
Affective_touch	0.38	0.38	-0.13	0.05	2.26	0.63
Sneeze	0.05	0.05	0.71	0.01	1.02	0.41
Cough	0.04	0.02	0.71	0.04	1.01	0.43
Wind	7.52e-03	-0.04	-5.18e-03	0.94	1.00	0.15
Burp	-0.01	0.27	0.17	0.47	1.91	0.48

The 4 latent factors (oblimin rotation) accounted for 44.25% of the total variance of the original data (MR1 = 17.71%, MR2 = 11.72%, MR3 = 7.65%, MR4 = 7.17%).

Sample 3

efa4 <- parameters::factor_analysis(data3, n = 4, rotation = "oblimin", sort = TRUE)
plot(efa4)

display(efa4)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR3	MR2	MR4	Complexity	Uniqueness
Cough	0.88	-0.08	-0.04	0.11	1.05	0.29
Sneeze	0.74	0.05	1.42e-03	-0.03	1.01	0.41
Vomit	0.55	0.12	0.04	-0.09	1.17	0.58
Burp	0.55	0.08	0.22	-0.13	1.47	0.48
Wind	0.53	0.10	0.21	-0.25	1.84	0.48
Temp	0.46	0.34	-0.06	0.09	1.97	0.52
Taste	0.29	0.18	0.10	0.12	2.37	0.74
Hungry	-8.85e-03	0.70	-0.01	-0.09	1.03	0.54
Urinate	0.08	0.66	-0.12	0.11	1.15	0.52
Thirsty	-0.08	0.64	-0.03	0.14	1.14	0.62
Defecate	0.10	0.63	-3.50e-05	-0.07	1.08	0.53
Breathing	-0.10	0.49	0.28	-0.07	1.76	0.66
Sex_arousal	0.17	0.44	0.17	-0.22	2.20	0.58
Pain	0.21	0.41	0.08	0.18	2.00	0.58
Heart	0.05	0.28	0.24	0.08	2.23	0.77
Muscles	0.24	0.27	0.20	0.19	3.66	0.61
Tickle	0.02	-0.05	0.84	-0.04	1.01	0.31
Itch	0.03	0.02	0.74	0.11	1.05	0.39
Affective_touch	0.08	0.23	0.36	0.12	2.08	0.67
Bruise	0.17	0.07	0.29	0.45	2.08	0.55
Blood_Sugar	4.11e-03	0.11	0.35	0.43	2.07	0.60

The 4 latent factors (oblimin rotation) accounted for 45.69% of the total variance of the original data (MR1 = 15.50%, MR3 = 15.46%, MR2 = 11.02%, MR4 = 3.71%).

Sample 4

efa5 <- parameters::factor_analysis(data4, n = 5, rotation = "oblimin", sort = TRUE)
plot(efa5)

display(efa5)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR3	MR5	MR4	MR2	Complexity	Uniqueness
Cough	0.67	-0.07	8.90e-03	0.06	0.03	1.04	0.54
Sneeze	0.63	-0.04	0.11	0.05	-0.08	1.11	0.54
Burp	0.58	0.18	-0.03	-0.09	0.03	1.25	0.59
Wind	0.54	0.08	-0.06	-0.03	0.11	1.17	0.67
Vomit	0.33	0.06	0.21	0.15	-0.19	3.02	0.70
Taste	0.24	0.19	-4.60e-03	0.06	0.18	3.00	0.77
Bruise	-0.03	0.63	0.05	0.05	-0.03	1.03	0.58
Affective_touch	0.07	0.53	0.03	-0.01	-0.04	1.05	0.67
Blood_Sugar	0.07	0.43	0.15	0.10	-0.18	1.82	0.67
Itch	0.11	0.41	-0.01	0.05	0.15	1.45	0.71
Muscles	0.07	0.40	-0.05	0.02	0.19	1.56	0.75
Urinate	-0.04	0.06	0.64	6.75e-03	0.06	1.04	0.55
Defecate	0.11	0.02	0.59	-0.03	0.02	1.08	0.57
Heart	-0.04	0.04	-0.03	0.66	-0.04	1.02	0.59
Breathing	0.09	-2.40e-03	-0.02	0.51	0.12	1.17	0.65
Sex_arousal	0.07	0.18	0.07	0.24	0.13	2.86	0.78
Thirsty	0.05	-0.05	0.24	0.03	0.49	1.51	0.60
Hungry	0.09	-0.10	0.15	0.21	0.38	2.25	0.67
Temp	0.09	0.17	7.51e-03	0.11	0.36	1.81	0.72
Pain	-3.56e-03	0.30	0.14	0.04	0.33	2.38	0.67

The 5 latent factors (oblimin rotation) accounted for 35.08% of the total variance of the original data (MR1 = 10.08%, MR3 = 8.61%, MR5 = 5.94%, MR4 = 5.33%, MR2 = 5.12%).

Sample 5

efa5 <- parameters::factor_analysis(data5, n = 5, rotation = "oblimin", sort = TRUE)
plot(efa5)

display(efa5)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR4	MR2	MR3	MR5	Complexity	Uniqueness
Burp	0.69	0.06	0.03	-0.10	0.03	1.06	0.50
Wind	0.62	0.01	0.03	-0.02	-0.08	1.04	0.64
Cough	0.54	-0.06	0.01	0.17	0.11	1.31	0.58
Sneeze	0.45	0.02	-0.07	0.21	0.17	1.78	0.61
Bruise	-0.04	0.60	-0.01	0.02	0.06	1.03	0.63
Affective_touch	0.08	0.53	-0.04	-0.08	0.06	1.12	0.67
Blood_Sugar	0.01	0.48	-0.09	0.08	0.19	1.45	0.68
Itch	0.18	0.38	0.08	0.12	-0.02	1.77	0.69
Muscles	0.09	0.38	0.19	3.00e-03	-0.07	1.68	0.77
Thirsty	0.03	-0.17	0.57	-0.01	0.22	1.50	0.55
Temp	2.56e-03	0.18	0.49	0.03	-0.08	1.34	0.71
Hungry	0.02	-0.02	0.48	0.16	0.09	1.32	0.65
Taste	0.19	0.22	0.32	0.05	-0.12	2.84	0.73
Pain	0.10	0.24	0.28	2.96e-03	0.03	2.26	0.77
Breathing	0.04	-0.03	0.03	0.67	-0.07	1.03	0.54
Heart	-0.08	0.05	5.68e-03	0.58	0.06	1.08	0.66
Sex_arousal	0.11	0.14	0.23	0.26	-0.06	3.08	0.76
Urinate	-0.04	0.11	0.11	-0.01	0.61	1.14	0.54
Defecate	0.14	0.03	0.03	-0.02	0.60	1.11	0.53
Vomit	0.15	0.15	-0.04	0.22	0.23	3.56	0.74

The 5 latent factors (oblimin rotation) accounted for 35.25% of the total variance of the original data (MR1 = 8.99%, MR4 = 8.07%, MR2 = 6.47%, MR3 = 5.97%, MR5 = 5.75%).

Sample 6

efa4 <- parameters::factor_analysis(data6, n = 4, rotation = "oblimin", sort = TRUE)
plot(efa4)

display(efa4)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR2	MR4	MR3	Complexity	Uniqueness
Hungry	0.69	0.22	-0.10	-0.15	1.36	0.47
Urinate	0.68	-0.13	-0.05	0.13	1.16	0.55
Thirsty	0.57	0.18	0.13	-0.08	1.34	0.47
Temp	0.56	-0.08	0.28	0.04	1.50	0.48
Pain	0.53	0.03	0.25	0.10	1.50	0.44
Defecate	0.49	-0.05	-0.05	0.31	1.75	0.60
Sex_arousal	0.46	0.07	-0.01	0.24	1.58	0.61
Affective_touch	0.43	0.18	0.18	-0.06	1.77	0.61
Muscles	0.36	0.02	0.26	0.26	2.71	0.51
Tickle	0.07	0.66	0.04	-2.12e-03	1.03	0.49
Itch	5.18e-04	0.65	-0.05	0.12	1.09	0.56
Bruise	0.02	0.59	0.11	-0.11	1.14	0.59
Burp	-0.02	0.44	0.04	0.32	1.86	0.61
Taste	0.13	0.33	0.17	0.32	2.82	0.50
Blood_Sugar	0.14	0.25	0.12	-0.09	2.44	0.85
Breathing	-0.03	0.03	0.81	2.25e-03	1.01	0.35
Heart	0.12	0.02	0.50	-8.27e-04	1.12	0.66
Vomit	0.19	0.14	0.27	0.16	3.06	0.66
Sneeze	0.03	0.28	0.12	0.45	1.84	0.54
Wind	0.18	0.14	0.15	0.42	1.88	0.54
Cough	0.14	0.21	0.22	0.25	3.53	0.62

The 4 latent factors (oblimin rotation) accounted for 44.31% of the total variance of the original data (MR1 = 16.65%, MR2 = 11.34%, MR4 = 9.45%, MR3 = 6.86%).

Sample 7

efa4 <- parameters::factor_analysis(data7, n = 4, rotation = "oblimin", sort = TRUE)
plot(efa4)

display(efa4)

Rotated loadings from Factor Analysis (oblimin-rotation)

Variable	MR1	MR3	MR4	MR2	Complexity	Uniqueness
Hungry	0.62	-0.04	0.15	0.04	1.13	0.51
Sex_arousal	0.58	0.06	0.01	0.03	1.03	0.60
Pain	0.55	0.04	0.09	-8.09e-03	1.06	0.61
Thirsty	0.55	0.09	0.07	0.03	1.10	0.58
Temp	0.36	0.11	0.19	0.09	1.85	0.65
Affective_touch	0.28	-0.05	0.13	0.26	2.51	0.75
Defecate	0.06	0.67	-7.30e-03	0.02	1.02	0.49
Sneeze	-0.18	0.53	0.29	0.06	1.85	0.59
Vomit	0.13	0.53	-0.05	0.09	1.20	0.61
Urinate	0.40	0.47	-0.12	-0.02	2.09	0.52
Burp	2.70e-03	0.46	0.17	0.08	1.34	0.64
Breathing	0.02	0.02	0.71	0.01	1.00	0.45
Heart	0.17	7.44e-03	0.59	0.01	1.16	0.51
Muscles	0.15	0.07	0.47	0.07	1.31	0.60
Itch	-5.06e-03	-0.02	-0.02	0.78	1.00	0.43
Tickle	-0.01	0.02	-0.01	0.73	1.00	0.47
Bruise	8.63e-03	0.04	0.06	0.38	1.08	0.81
Wind	0.10	0.08	0.10	0.26	1.82	0.82

The 4 latent factors (oblimin rotation) accounted for 40.83% of the total variance of the original data (MR1 = 12.34%, MR3 = 10.09%, MR4 = 9.29%, MR2 = 9.12%).

CFA

Structure

All data

# Ambiguous: Temp + Vomit + Affective_touch + Sex_arousal + Taste
# Discard: Tickle + Blood_Sugar

m1 <- "
Interoception =~ Hungry + Thirsty + Urinate + Defecate + Itch + Bruise + Muscles + Pain + Breathing + Heart + Cough + Sneeze + Wind + Burp
"

m4 <- "
Sustenance =~ Hungry + Thirsty + Urinate + Defecate + Muscles + Pain
Skin =~ Itch + Bruise
HeartBreathing =~ Breathing + Heart
Expulsion =~ Cough + Sneeze + Wind + Burp
"


m5 <- "
Sustenance =~ Hungry + Thirsty + Pain
UrinateDefecate =~ Urinate + Defecate
Skin =~ Itch + Bruise + Muscles
HeartBreathing =~ Breathing + Heart
Expulsion =~ Cough + Sneeze + Wind + Burp
"


m6 <- "
HungryThirsty =~ Hungry + Thirsty
UrinateDefecate =~ Urinate + Defecate
ItchBruise =~ Itch + Bruise
MusclesPain =~ Muscles + Pain
HeartBreathing =~ Breathing + Heart
Expulsion =~ Cough + Sneeze + Wind + Burp
"

m7 <- "
HungryThirsty =~ Hungry + Thirsty
UrinateDefecate =~ Urinate + Defecate
ItchBruise =~ Itch + Bruise
MusclesPain =~ Muscles + Pain
HeartBreathing =~ Breathing + Heart
CoughSneeze =~ Cough + Sneeze
WindBurp =~ Wind + Burp
"

m1_all <- lavaan::cfa(m1, data = data_all)
m4_all <- lavaan::cfa(m4, data = data_all)
m5_all <- lavaan::cfa(m5, data = data_all)
m6_all <- lavaan::cfa(m6, data = data_all)
m7_all <- lavaan::cfa(m7, data = data_all)

anova(m1_all, m4_all, m5_all, m6_all, m7_all) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
m7_all	-14218.06	-13919.70	239.58	56
m6_all	-14036.69	-13774.87	432.95	62	< .001
m5_all	-13745.07	-13513.69	734.56	67	< .001
m4_all	-13389.39	-13182.37	1098.24	71	< .001
m1_all	-12372.54	-12202.05	2127.10	77	< .001

Anova Table (Type 1 tests)

parameters::parameters(m7_all, standardize=TRUE) |>
  display()

# Loading

Link	Coefficient	SE	95% CI	z	p
HungryThirsty =~ Hungry	0.66	0.02	(0.63, 0.70)	40.45	< .001
HungryThirsty =~ Thirsty	0.69	0.02	(0.66, 0.73)	42.07	< .001
UrinateDefecate =~ Urinate	0.69	0.02	(0.66, 0.72)	45.07	< .001
UrinateDefecate =~ Defecate	0.71	0.02	(0.68, 0.74)	46.26	< .001
ItchBruise =~ Itch	0.59	0.02	(0.55, 0.62)	32.53	< .001
ItchBruise =~ Bruise	0.56	0.02	(0.52, 0.59)	30.99	< .001
MusclesPain =~ Muscles	0.58	0.02	(0.55, 0.61)	36.37	< .001
MusclesPain =~ Pain	0.68	0.02	(0.64, 0.71)	42.39	< .001
HeartBreathing =~ Breathing	0.70	0.02	(0.66, 0.74)	34.75	< .001
HeartBreathing =~ Heart	0.56	0.02	(0.52, 0.60)	29.78	< .001
CoughSneeze =~ Cough	0.73	0.01	(0.70, 0.75)	55.74	< .001
CoughSneeze =~ Sneeze	0.71	0.01	(0.69, 0.74)	54.21	< .001
WindBurp =~ Wind	0.67	0.01	(0.64, 0.69)	46.98	< .001
WindBurp =~ Burp	0.74	0.01	(0.71, 0.77)	53.17	< .001

# Correlation

Link	Coefficient	SE	95% CI	z	p
HungryThirsty ~~ UrinateDefecate	0.59	0.02	(0.55, 0.64)	26.16	< .001
HungryThirsty ~~ ItchBruise	0.49	0.03	(0.43, 0.54)	16.82	< .001
HungryThirsty ~~ MusclesPain	0.61	0.02	(0.56, 0.66)	24.64	< .001
HungryThirsty ~~ HeartBreathing	0.51	0.03	(0.46, 0.57)	19.82	< .001
HungryThirsty ~~ CoughSneeze	0.47	0.02	(0.43, 0.52)	20.06	< .001
HungryThirsty ~~ WindBurp	0.44	0.02	(0.39, 0.49)	17.95	< .001
UrinateDefecate ~~ ItchBruise	0.45	0.03	(0.39, 0.51)	15.70	< .001
UrinateDefecate ~~ MusclesPain	0.63	0.02	(0.58, 0.67)	26.14	< .001
UrinateDefecate ~~ HeartBreathing	0.48	0.03	(0.43, 0.53)	18.57	< .001
UrinateDefecate ~~ CoughSneeze	0.55	0.02	(0.51, 0.59)	25.07	< .001
UrinateDefecate ~~ WindBurp	0.50	0.02	(0.45, 0.54)	21.56	< .001
ItchBruise ~~ MusclesPain	0.82	0.03	(0.77, 0.88)	29.25	< .001
ItchBruise ~~ HeartBreathing	0.52	0.03	(0.46, 0.58)	17.22	< .001
ItchBruise ~~ CoughSneeze	0.64	0.03	(0.59, 0.69)	24.48	< .001
ItchBruise ~~ WindBurp	0.67	0.03	(0.62, 0.73)	25.64	< .001
MusclesPain ~~ HeartBreathing	0.56	0.03	(0.50, 0.61)	20.38	< .001
MusclesPain ~~ CoughSneeze	0.57	0.02	(0.52, 0.62)	23.69	< .001
MusclesPain ~~ WindBurp	0.59	0.02	(0.54, 0.63)	24.12	< .001
HeartBreathing ~~ CoughSneeze	0.51	0.02	(0.46, 0.55)	20.42	< .001
HeartBreathing ~~ WindBurp	0.42	0.03	(0.37, 0.47)	16.21	< .001
CoughSneeze ~~ WindBurp	0.76	0.02	(0.73, 0.80)	41.66	< .001

Sample 1

anova(
  lavaan::cfa(m1, data = data1),
  lavaan::cfa(m4, data = data1),
  lavaan::cfa(m5, data = data1),
  lavaan::cfa(m6, data = data1),
  lavaan::cfa(m7, data = data1)
  ) |>
  parameters::parameters() |>
  display()

Warning in lav_object_post_check(object): lavaan WARNING: covariance matrix of latent variables
                is not positive definite;
                use lavInspect(fit, "cov.lv") to investigate.

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data1)	-1950.88	-1749.41	88.22	56
lavaan::cfa(m6, data = data1)	-1924.80	-1748.01	126.29	62	< .001
lavaan::cfa(m5, data = data1)	-1859.71	-1703.48	201.38	67	< .001
lavaan::cfa(m4, data = data1)	-1834.28	-1694.49	234.82	71	< .001
lavaan::cfa(m1, data = data1)	-1660.86	-1545.74	420.24	77	< .001

Anova Table (Type 1 tests)

Sample 2

anova(
  lavaan::cfa(m1, data = data2),
  lavaan::cfa(m4, data = data2),
  lavaan::cfa(m5, data = data2),
  lavaan::cfa(m6, data = data2),
  lavaan::cfa(m7, data = data2)
  ) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data2)	-1898.38	-1705.96	66.50	56
lavaan::cfa(m6, data = data2)	-1822.70	-1653.84	154.19	62	< .001
lavaan::cfa(m5, data = data2)	-1800.91	-1651.68	185.98	67	< .001
lavaan::cfa(m4, data = data2)	-1754.66	-1621.14	240.22	71	< .001
lavaan::cfa(m1, data = data2)	-1632.46	-1522.51	374.42	77	< .001

Anova Table (Type 1 tests)

Sample 3

anova(
  lavaan::cfa(m1, data = data3),
  lavaan::cfa(m4, data = data3),
  lavaan::cfa(m5, data = data3),
  lavaan::cfa(m6, data = data3),
  lavaan::cfa(m7, data = data3)
  ) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data3)	-3587.97	-3380.01	114.03	56
lavaan::cfa(m6, data = data3)	-3507.56	-3325.06	206.45	62	< .001
lavaan::cfa(m5, data = data3)	-3452.59	-3291.31	271.41	67	< .001
lavaan::cfa(m4, data = data3)	-3418.95	-3274.65	313.05	71	< .001
lavaan::cfa(m1, data = data3)	-3143.19	-3024.35	600.82	77	< .001

Anova Table (Type 1 tests)

Sample 4

anova(
  lavaan::cfa(m1, data = data4),
  lavaan::cfa(m4, data = data4),
  lavaan::cfa(m5, data = data4),
  lavaan::cfa(m6, data = data4),
  lavaan::cfa(m7, data = data4)
  ) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data4)	-5136.94	-4888.93	115.89	56
lavaan::cfa(m6, data = data4)	-5103.03	-4885.39	161.81	62	< .001
lavaan::cfa(m5, data = data4)	-5041.07	-4848.74	233.76	67	< .001
lavaan::cfa(m4, data = data4)	-4930.21	-4758.13	352.62	71	< .001
lavaan::cfa(m1, data = data4)	-4591.03	-4449.31	703.81	77	< .001

Anova Table (Type 1 tests)

Sample 5

anova(
  lavaan::cfa(m1, data = data5),
  lavaan::cfa(m4, data = data5),
  lavaan::cfa(m5, data = data5),
  lavaan::cfa(m6, data = data5),
  lavaan::cfa(m7, data = data5)
  ) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data5)	-1737.85	-1531.33	83.61	56
lavaan::cfa(m6, data = data5)	-1711.33	-1530.10	122.14	62	< .001
lavaan::cfa(m5, data = data5)	-1683.25	-1523.10	160.21	67	< .001
lavaan::cfa(m4, data = data5)	-1629.89	-1486.59	221.58	71	< .001
lavaan::cfa(m1, data = data5)	-1485.55	-1367.54	377.91	77	< .001

Anova Table (Type 1 tests)

Sample 6

anova(
  lavaan::cfa(m1, data = data6),
  lavaan::cfa(m4, data = data6),
  lavaan::cfa(m5, data = data6),
  lavaan::cfa(m6, data = data6),
  lavaan::cfa(m7, data = data6)
  ) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data6)	-1437.92	-1265.17	91.47	56
lavaan::cfa(m6, data = data6)	-1444.97	-1293.38	96.42	62	0.551
lavaan::cfa(m5, data = data6)	-1412.50	-1278.53	138.89	67	< .001
lavaan::cfa(m4, data = data6)	-1419.48	-1299.62	139.91	71	0.907
lavaan::cfa(m1, data = data6)	-1342.29	-1243.58	229.10	77	< .001

Anova Table (Type 1 tests)

Sample 7

anova(
  lavaan::cfa(str_remove(m5, fixed("Cough + ")), data = data7),
  lavaan::cfa(str_remove(m1, fixed("Cough + ")), data = data7),
  lavaan::cfa(str_remove(m4, fixed("Cough + ")), data = data7),
  lavaan::cfa(str_remove(m6, fixed("Cough + ")), data = data7),
  lavaan::cfa(str_remove(m7, fixed("\nCoughSneeze =~ Cough + Sneeze")), data = data7)
  ) |>
  parameters::parameters() |>
  display()

Warning in lav_object_post_check(object): lavaan WARNING: covariance matrix of latent variables
                is not positive definite;
                use lavInspect(fit, "cov.lv") to investigate.

Warning in lavTestLRT(object = object, ..., model.names = NAMES): lavaan
WARNING: some models are based on a different set of observed variables

Parameter	AIC	BIC	Chi2	df	p
c(“lavaan::cfa(str_remove(m7, fixed("\nCoughSneeze =~ Cough + Sneeze")),”, ” data = data7)“)	-1606.97	-1443.78	82.59	39
lavaan::cfa(str_remove(m6, fixed(“Cough +”)), data = data7)	-1825.70	-1654.15	109.52	50	0.005
lavaan::cfa(str_remove(m5, fixed(“Cough +”)), data = data7)	-1821.75	-1671.12	123.48	55	0.016
lavaan::cfa(str_remove(m4, fixed(“Cough +”)), data = data7)	-1752.45	-1618.56	200.78	59	< .001
lavaan::cfa(str_remove(m1, fixed(“Cough +”)), data = data7)	-1665.13	-1556.35	300.09	65	< .001

Anova Table (Type 1 tests)

Higher Order Factors

All data

m7h1 <- paste0(m7, "
              Interoception =~ HungryThirsty + UrinateDefecate + ItchBruise + MusclesPain + HeartBreathing + CoughSneeze + WindBurp
              ")

m7h3 <- paste0(m7, "
              Valenced =~ MusclesPain + HeartBreathing + ItchBruise
              Expulsion =~ CoughSneeze + WindBurp
              Homeostasis =~ HungryThirsty + UrinateDefecate
              ")

m7h1_all <- lavaan::cfa(m7h1, data = data_all)
m7h3_all <- lavaan::cfa(m7h3, data = data_all)
anova(m7_all, m7h1_all, m7h3_all) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
m7_all	-14218.06	-13919.70	239.58	56
m7h3_all	-14104.40	-13873.02	375.24	67	< .001
m7h1_all	-13903.22	-13690.11	582.42	70	< .001

Anova Table (Type 1 tests)

Sample 1

anova(
  lavaan::cfa(m7, data = data1), 
  lavaan::cfa(m7h1, data = data1), 
  lavaan::cfa(m7h3, data = data1)) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data1)	-1950.88	-1749.41	88.22	56
lavaan::cfa(m7h3, data = data1)	-1903.59	-1747.35	157.51	67	< .001
lavaan::cfa(m7h1, data = data1)	-1873.38	-1729.48	193.71	70	< .001

Anova Table (Type 1 tests)

Sample 2

anova(
  lavaan::cfa(m7, data = data2), 
  lavaan::cfa(m7h1, data = data2), 
  lavaan::cfa(m7h3, data = data2)) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data2)	-1898.38	-1705.96	66.50	56
lavaan::cfa(m7h3, data = data2)	-1855.12	-1705.90	131.76	67	< .001
lavaan::cfa(m7h1, data = data2)	-1841.29	-1703.85	151.59	70	< .001

Anova Table (Type 1 tests)

Sample 3

anova(
  lavaan::cfa(m7, data = data3), 
  lavaan::cfa(m7h1, data = data3), 
  lavaan::cfa(m7h3, data = data3)) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data3)	-3587.97	-3380.01	114.03	56
lavaan::cfa(m7h3, data = data3)	-3580.71	-3419.43	143.30	67	0.002
lavaan::cfa(m7h1, data = data3)	-3496.70	-3348.16	233.30	70	< .001

Anova Table (Type 1 tests)

Sample 4

anova(
  lavaan::cfa(m7, data = data4), 
  lavaan::cfa(m7h1, data = data4), 
  lavaan::cfa(m7h3, data = data4)) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data4)	-5136.94	-4888.93	115.89	56
lavaan::cfa(m7h3, data = data4)	-5103.39	-4911.06	171.44	67	< .001
lavaan::cfa(m7h1, data = data4)	-5007.06	-4829.92	273.77	70	< .001

Anova Table (Type 1 tests)

Sample 5

anova(
  lavaan::cfa(m7, data = data5), 
  lavaan::cfa(m7h1, data = data5), 
  lavaan::cfa(m7h3, data = data5)) |>
  parameters::parameters() |>
  display()

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data5)	-1737.85	-1531.33	83.61	56
lavaan::cfa(m7h3, data = data5)	-1727.18	-1567.03	116.28	67	< .001
lavaan::cfa(m7h1, data = data5)	-1709.86	-1562.35	139.61	70	< .001

Anova Table (Type 1 tests)

Sample 6

anova(
  lavaan::cfa(m7, data = data6), 
  lavaan::cfa(m7h1, data = data6), 
  lavaan::cfa(m7h3, data = data6)) |>
  parameters::parameters() |>
  display()

Warning in lav_object_post_check(object): lavaan WARNING: covariance matrix of latent variables
                is not positive definite;
                use lavInspect(fit, "cov.lv") to investigate.

Parameter	AIC	BIC	Chi2	df	p
lavaan::cfa(m7, data = data6)	-1437.92	-1265.17	91.47	56
lavaan::cfa(m7h3, data = data6)	-1438.54	-1304.58	112.85	67	0.030
lavaan::cfa(m7h1, data = data6)	-1430.27	-1306.88	127.12	70	0.003

Anova Table (Type 1 tests)

Sample 7

None of the model can be properly estimated (variances are negative).

NA

[1] NA

# anova(
#   lavaan::cfa(str_remove(m7, fixed("\nCoughSneeze =~ Cough + Sneeze")), data = data7), 
#   lavaan::cfa(str_remove(str_remove(m7h1, fixed("\nCoughSneeze =~ Cough + Sneeze")), fixed("CoughSneeze +")), data = data7), 
#   lavaan::cfa(str_remove(str_remove(m7h3, fixed("\nCoughSneeze =~ Cough + Sneeze")), fixed("CoughSneeze +")), data = data7)) |>
#   parameters::parameters() |>
#   display()

Model Performance

rbind(
  performance::performance(m7_all, metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "All"),
  performance::performance(update(m7_all, data = data1), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 1"),
  performance::performance(update(m7_all, data = data2), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 2"),
  performance::performance(update(m7_all, data = data3), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 3"),
  performance::performance(update(m7_all, data = data4), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 4"),
  performance::performance(update(m7_all, data = data5), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 5"),
  performance::performance(update(m7_all, data = data6), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 6"),
  performance::performance(lavaan::cfa(str_remove(m7, fixed("\nCoughSneeze =~ Cough + Sneeze")), data = data7), metrics = c("Chi2", "RMSEA", "CFI", "SRMR")) |>
    mutate(Sample = "Sample 7")
) |>
  display()

Chi2	RMSEA	CFI	SRMR	Sample
239.58	0.03	0.98	0.02	All
88.22	0.04	0.98	0.03	Sample 1
66.50	0.02	0.99	0.03	Sample 2
114.03	0.04	0.98	0.03	Sample 3
115.89	0.03	0.98	0.02	Sample 4
83.61	0.03	0.98	0.03	Sample 5
91.47	0.05	0.97	0.04	Sample 6
82.59	0.05	0.97	0.03	Sample 7

No evidence for higher order factors.

Scores

Note the usage of descriptive factor names relating directly to the items to avoid abstraction.

scores <- list(
  sample1 = cbind(ega_scores_1, as.data.frame(predict(lavaan::cfa(m7, data = data1)))) |>
    mutate(Sample = "Sample 1"),
  sample2 = cbind(ega_scores_2, as.data.frame(predict(lavaan::cfa(m7, data = data2)))) |>
    mutate(Sample = "Sample 2"),
  sample3 = cbind(ega_scores_3, as.data.frame(predict(lavaan::cfa(m7, data = data3)))) |>
    mutate(Sample = "Sample 3"),
  sample4 = cbind(ega_scores_4, as.data.frame(predict(lavaan::cfa(m7, data = data4)))) |>
    mutate(Sample = "Sample 4"),
  sample5 = cbind(ega_scores_5, as.data.frame(predict(lavaan::cfa(m7, data = data5)))) |>
    mutate(Sample = "Sample 5"),
  sample6 = cbind(ega_scores_6, as.data.frame(predict(lavaan::cfa(m7, data = data6)))) |>
    mutate(Sample = "Sample 6"),
  sample7 = cbind(ega_scores_7, as.data.frame(predict(lavaan::cfa(str_remove(m7, fixed("\nCoughSneeze =~ Cough + Sneeze")), data = data7)))) |>
    mutate(Sample = "Sample 7", CoughSneeze = NA)
)
  
save(scores, file = "../data/scores.RData")

Study 2 - Correlates

load("../data/scores.RData")

df1 <- cbind(df1, scores$sample1)
df2 <- cbind(df2, scores$sample2)
df3 <- cbind(df3, scores$sample3)
df4 <- cbind(df4, scores$sample4)
df5 <- cbind(df5, scores$sample5)
df6 <- cbind(df6, scores$sample6)
df7 <- cbind(df7, scores$sample7)

vars_intero <- names(select(scores$sample1, -Sample))
vars_cfa <- names(select(scores$sample1, -Sample, -starts_with("EGA")))
vars_cfa7 <- vars_cfa[vars_cfa != "CoughSneeze"]

# names(df1)
# names(df2)
# names(df3)  # TAS, BPQ
# names(df4)  # BPQ
# names(df5)  # PHQ, TAS
# names(df6)  # BDI, STAI
# names(df7)  # MAIA, PHQ4, PI, LIE, IPIP, PID, SPQ

EGA and CFA

df_scores <- rbind(
  select(df1, all_of(vars_intero)),
  select(df2, all_of(vars_intero)),
  select(df3, all_of(vars_intero)),
  select(df4, all_of(vars_intero)),
  select(df5, all_of(vars_intero)),
  select(df6, all_of(vars_intero)),
  select(df7, all_of(vars_intero))
)

correlation::correlation(df_scores, redundant = TRUE) |>
    correlation::cor_sort() |>
    correlation::cor_lower() |>
    mutate(val = paste0(insight::format_value(r), format_p(p, stars_only = TRUE))) |>
    mutate(Parameter2 = fct_rev(Parameter2)) |>
    ggplot(aes(x = Parameter1, y = Parameter2)) +
    geom_tile(aes(fill = r), color = "white") +
    geom_text(aes(label = val), size = 3) +
    labs(title = "Correlation Matrix") +
    scale_fill_metro_c(limit = c(0, 1), guide = guide_colourbar(ticks = FALSE)) +
    theme_minimal() +
    theme(
      legend.title = element_blank(),
      axis.title.x = element_blank(),
      axis.title.y = element_blank(),
      axis.text.x = element_text(angle = 45, hjust = 1)
    )

Demographics

make_diff <- function(df, vars, component=NA) {
  dat <- data.frame()
  for (resp in vars) {
    f <- as.formula(paste0(resp, "~ Gender"))
    d <- filter(df, Gender %in% c("Male", "Female")) |> 
      mutate(Gender = fct_relevel(Gender, "Male", "Female"))
    dat <- effectsize::cohens_d(f, data=d) |> 
      mutate(var_Correlate = "Gender", var_Interoception = resp, 
             p=t.test(f, data=d)$p.val) |> 
      as.data.frame() |> 
      rbind(dat)
  }
    
  dat$r <- effectsize::d_to_r(dat$Cohens_d, n1=sum(df$Gender=="Male"), n2=sum(df$Gender=="Female"))
  dat$CI_low <- effectsize::d_to_r(dat$CI_low, n1=sum(df$Gender=="Male"), n2=sum(df$Gender=="Female"))
  dat$CI_high <- effectsize::d_to_r(dat$CI_high, n1=sum(df$Gender=="Male"), n2=sum(df$Gender=="Female"))
  dat$Sample <- unique(df$Sample)
  dat$Component <- component
  select(dat, -Cohens_d)
}

data <- rbind(
  make_diff(df1, vars_cfa, component="Demographic"),
  make_diff(df2, vars_cfa, component="Demographic"),
  make_diff(df3, vars_cfa, component="Demographic"),
  make_diff(df4, vars_cfa, component="Demographic"),
  make_diff(df5, vars_cfa, component="Demographic"),
  make_diff(df6, vars_cfa, component="Demographic"),
  make_diff(df7, vars_cfa7, component="Demographic")
) |> 
  mutate(var_Correlate = "Gender (Male vs. Female)")

make_cor <- function(df, vars, with = "Age", component=NA) {
  dat <- data.frame() 
  for (v in vars) {
    dat <- correlation::cor_test(df, v, with) |> 
      rename(var_Interoception = Parameter1, var_Correlate = Parameter2) |> 
      as.data.frame() |> 
      rbind(dat)
  }
    
  dat$Sample <- unique(df$Sample)
  dat$Component <- component
  select(dat, -df_error, -t, -Method, -n_Obs)
}

data <- rbind(
  data,
  make_cor(df1, vars_cfa, with = "Age", component="Demographic"),
  make_cor(df2, vars_cfa, with = "Age", component="Demographic"),
  make_cor(df3, vars_cfa, with = "Age", component="Demographic"),
  make_cor(df4, vars_cfa, with = "Age", component="Demographic"),
  make_cor(df5, vars_cfa, with = "Age", component="Demographic"),
  make_cor(df6, vars_cfa, with = "Age", component="Demographic"),
  make_cor(df7, vars_cfa7, with = "Age", component="Demographic")
) 

Body

# Reverse (TODO: Ana confirm why)
df3$BPQ_BodyAwareness <- -1*normalize(df3$BPQ_BodyAwareness)
df3$BPQ_AutonomicReactivity <- -1*normalize(df3$BPQ_AutonomicReactivity)

data <- rbind(
  data,
  make_cor(df3, vars_cfa, with = "BPQ_BodyAwareness", component="Body") |> 
    mutate(var_Correlate = "Body Awareness"),
  make_cor(df3, vars_cfa, with = "BPQ_AutonomicReactivity", component="Body") |> 
    mutate(var_Correlate = "Autonomic Reactivity"),
  make_cor(df4, vars_cfa, with = "BPQ_A", component="Body") |> 
    mutate(var_Correlate = "Body Awareness")
)

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "MAIA_AttentionRegulation", component="MAIA")  |> 
    mutate(var_Correlate = "Attention Regulation"),
  make_cor(df7, vars_cfa7, with = "MAIA_BodyListening", component="MAIA") |> 
    mutate(var_Correlate = "Body Listening"),
  make_cor(df7, vars_cfa7, with = "MAIA_EmotionalAwareness", component="MAIA") |> 
    mutate(var_Correlate = "Emotional Awareness"),
  make_cor(df7, vars_cfa7, with = "MAIA_NotDistracting", component="MAIA") |> 
    mutate(var_Correlate = "Not Distracting"),
  make_cor(df7, vars_cfa7, with = "MAIA_Noticing", component="MAIA") |> 
    mutate(var_Correlate = "Noticing"),
  make_cor(df7, vars_cfa7, with = "MAIA_NotWorrying", component="MAIA") |> 
    mutate(var_Correlate = "Not Worrying"),
  make_cor(df7, vars_cfa7, with = "MAIA_SelfRegulation", component="MAIA") |> 
    mutate(var_Correlate = "Self Regulation"),
  make_cor(df7, vars_cfa7, with = "MAIA_Trusting", component="MAIA") |> 
    mutate(var_Correlate = "Trusting")
)

# names(df3)[str_detect(names(df3), "TAS")]
# names(df5)[str_starts(names(df5), "TAS")]

data <- rbind(
  data,
  make_cor(df3, vars_cfa, with = "TAS_IdentifyingFeelings", component="Body") |> 
    mutate(var_Correlate = "Difficulties Identifying Feelings"),
  make_cor(df3, vars_cfa, with = "TAS_DescribingFeelings", component="Body") |> 
    mutate(var_Correlate = "Difficulties Describing Feelings"),
  make_cor(df3, vars_cfa, with = "TAS_ExternallyOrientedThinking", component="Body") |> 
    mutate(var_Correlate = "External Thinking"),
  make_cor(df5, vars_cfa, with = "TAS_DIF", component="Body") |> 
    mutate(var_Correlate = "Difficulties Identifying Feelings"),
  make_cor(df5, vars_cfa, with = "TAS_DDF", component="Body") |> 
    mutate(var_Correlate = "Difficulties Describing Feelings")
)

Mood

data <- rbind(
  data,
  make_cor(df5, vars_cfa, with = "PHQ15_sum", component="Mood") |> 
    mutate(var_Correlate = "Somatic Concerns (PHQ-15)"),
  make_cor(df5, vars_cfa, with = "PHQ9_sum", component="Mood")|> 
    mutate(var_Correlate = "Depression (PHQ-9)"),
  make_cor(df7, vars_cfa7, with = "PHQ4_Depression", component="Mood") |> 
    mutate(var_Correlate = "Depression (PHQ-4)"),
  make_cor(df7, vars_cfa7, with = "PHQ4_Anxiety", component="Mood") |> 
    mutate(var_Correlate = "Anxiety (PHQ-4)")
)

data <- rbind(
  data,
  make_cor(df6, vars_cfa, with = "BDI2_Total", component="Mood") |> 
    mutate(var_Correlate = "Depression (BDI-II)")
)

data <- rbind(
  data,
  make_cor(df6, vars_cfa, with = "STAI5_General", component="Mood") |> 
    mutate(var_Correlate = "Anxiety (STAI-5)")
)

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "BPD", component="Mood") |> 
    mutate(var_Correlate = "Boderline Personality")
)

Psychopathology

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "PID5_Antagonism", component="Maladaptive") |> 
    mutate(var_Correlate = "Antagonism"),
  make_cor(df7, vars_cfa7, with = "PID5_Detachment", component="Maladaptive") |> 
    mutate(var_Correlate = "Detachment"),
  make_cor(df7, vars_cfa7, with = "PID5_Disinhibition", component="Maladaptive") |> 
    mutate(var_Correlate = "Disinhibition"),
  make_cor(df7, vars_cfa7, with = "PID5_NegativeAffect", component="Maladaptive") |> 
    mutate(var_Correlate = "Negative Affect"),
  make_cor(df7, vars_cfa7, with = "PID5_Psychoticism", component="Maladaptive") |> 
    mutate(var_Correlate = "Psychoticism")
)

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "SPQ_ConstrictedAffect", component="Schizotypic") |> 
    mutate(var_Correlate = "Constricted Affect"),
  make_cor(df7, vars_cfa7, with = "SPQ_Eccentric", component="Schizotypic") |> 
    mutate(var_Correlate = "Eccentric"),
  make_cor(df7, vars_cfa7, with = "SPQ_MagicalThinking", component="Schizotypic") |> 
    mutate(var_Correlate = "Magical Thinking"),
  make_cor(df7, vars_cfa7, with = "SPQ_NoCloseFriends", component="Schizotypic") |> 
    mutate(var_Correlate = "No Close Friends"),
  make_cor(df7, vars_cfa7, with = "SPQ_OddSpeech", component="Schizotypic") |> 
    mutate(var_Correlate = "Odd Speech"),
  make_cor(df7, vars_cfa7, with = "SPQ_Reference", component="Schizotypic") |> 
    mutate(var_Correlate = "Reference"),
  make_cor(df7, vars_cfa7, with = "SPQ_SocialAnxiety", component="Schizotypic") |> 
    mutate(var_Correlate = "Social Anxiety"),
  make_cor(df7, vars_cfa7, with = "SPQ_Suspiciousness", component="Schizotypic") |> 
    mutate(var_Correlate = "Suspiciousness"),
  make_cor(df7, vars_cfa7, with = "SPQ_UnusualPerceptions", component="Schizotypic") |> 
    mutate(var_Correlate = "Unusual Perceptions")
)

#  names(select(df7, starts_with("ASQ")))
data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "ASQ_Imagination", component="Autistic") |> 
    mutate(var_Correlate = "Imagination"),
  make_cor(df7, vars_cfa7, with = "ASQ_LackSocialSkills", component="Autistic") |>
    mutate(var_Correlate = "Lack of Social Skills"),
  make_cor(df7, vars_cfa7, with = "ASQ_LowAttentionalSwitching", component="Autistic") |>
    mutate(var_Correlate = "Low Attentional Switching"),
  make_cor(df7, vars_cfa7, with = "ASQ_Patterns", component="Autistic") |>
    mutate(var_Correlate = "Patterns and Numbers"),
  make_cor(df7, vars_cfa7, with = "ASQ_Routine", component="Autistic") |>
    mutate(var_Correlate = "Routines")
)

Personality

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "IPIP6_Agreeableness", component="Personality") |> 
    mutate(var_Correlate = "Agreeableness"),
  make_cor(df7, vars_cfa7, with = "IPIP6_Conscientiousness", component="Personality") |>
    mutate(var_Correlate = "Conscientiousness"),
  make_cor(df7, vars_cfa7, with = "IPIP6_Neuroticism", component="Personality") |>
    mutate(var_Correlate = "Neuroticism"),
  make_cor(df7, vars_cfa7, with = "IPIP6_Openness", component="Personality") |>
    mutate(var_Correlate = "Openness"),
  make_cor(df7, vars_cfa7, with = "IPIP6_Extraversion", component="Personality") |>
    mutate(var_Correlate = "Extraversion"),
  make_cor(df7, vars_cfa7, with = "IPIP6_HonestyHumility", component="Personality") |>
    mutate(var_Correlate = "Honesty-Humility")
)

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "PI_Alive", component="World Beliefs") |> 
    mutate(var_Correlate = "Alive"),
  make_cor(df7, vars_cfa7, with = "PI_Changing", component="World Beliefs") |>
    mutate(var_Correlate = "Changing"),
  make_cor(df7, vars_cfa7, with = "PI_Enticing", component="World Beliefs") |>
    mutate(var_Correlate = "Enticing"),
  make_cor(df7, vars_cfa7, with = "PI_Good", component="World Beliefs") |>
    mutate(var_Correlate = "Good"),
  make_cor(df7, vars_cfa7, with = "PI_Hierarchical", component="World Beliefs") |>
    mutate(var_Correlate = "Hierarchical"),
  make_cor(df7, vars_cfa7, with = "PI_Safe", component="World Beliefs") |>
    mutate(var_Correlate = "Safe")
)

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "GCBS_Extraterrestrial", component="Conspiracy Beliefs") |> 
    mutate(var_Correlate = "Extraterrestrial"),
  make_cor(df7, vars_cfa7, with = "GCBS_GlobalConspiracies", component="Conspiracy Beliefs") |>
    mutate(var_Correlate = "Global Conspiracies"),
  make_cor(df7, vars_cfa7, with = "GCBS_GovernmentMalfeasance", component="Conspiracy Beliefs") |>
    mutate(var_Correlate = "Government Malfeasance"),
  make_cor(df7, vars_cfa7, with = "GCBS_InformationControl", component="Conspiracy Beliefs") |>
    mutate(var_Correlate = "Information Control"),
  make_cor(df7, vars_cfa7, with = "GCBS_PersonalWellbeing", component="Conspiracy Beliefs") |>
    mutate(var_Correlate = "Personal Wellbeing")
)

data <- rbind(
  data,
  make_cor(df7, vars_cfa7, with = "LIE_Ability", component="Lying Profile") |> 
    mutate(var_Correlate = "Ability"),
  make_cor(df7, vars_cfa7, with = "LIE_Frequency", component="Lying Profile") |>
    mutate(var_Correlate = "Frequency"),
  make_cor(df7, vars_cfa7, with = "LIE_Contextuality", component="Lying Profile") |>
    mutate(var_Correlate = "Contextuality"),
  make_cor(df7, vars_cfa7, with = "LIE_Negativity", component="Lying Profile") |>
    mutate(var_Correlate = "Negativity")
)

Summary

p <- data |>
  mutate(sig = ifelse(p < .01, "p < .01", ifelse(p < .001, "p < .001", "N.S.")),
         dir = sign(r),
         var_Interoception = str_replace(var_Interoception, "CoughSneeze", "Cough/Sneeze"),
         var_Interoception = str_replace(var_Interoception, "HeartBreathing", "Heart/Breathing"),
         var_Interoception = str_replace(var_Interoception, "HungryThirsty", "Hungry/Thirsty"),
         var_Interoception = str_replace(var_Interoception, "ItchBruise", "Itch/Bruise"),
         var_Interoception = str_replace(var_Interoception, "MusclesPain", "Muscles/Pain"),
         var_Interoception = str_replace(var_Interoception, "UrinateDefecate", "Urinate/Defecate"),
         var_Interoception = str_replace(var_Interoception, "WindBurp", "Wind/Burp"),
         var_Interoception = str_replace(var_Interoception, "Expulsion_Sudden", "Expulsion - Sudden"),
         var_Interoception = fct_relevel(var_Interoception, "Heart/Breathing", "Muscles/Pain", "Itch/Bruise", "Hungry/Thirsty", "Urinate/Defecate", "Wind/Burp", "Cough/Sneeze"),
         Component = fct_relevel(
           Component, "Demographic", "MAIA", "Body", "Mood", "Maladaptive",
           "Schizotypic", "Autistic", "Personality", "Conspiracy Beliefs", "World Beliefs", "Lying Profile")) |>
  ggplot(aes(x = var_Correlate)) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  # geom_bar(aes(fill=Sample, alpha=sig, y=r, group=Sample), stat="identity", position = position_dodge2(width = 0.8)) +
  geom_linerange(aes(color=dir, alpha=sig, ymin=0, ymax=r, group=Sample), position = position_dodge2(width = 0.7), linewidth=2.5) +
  geom_linerange(aes(ymin = CI_low, ymax = CI_high, color = dir, alpha=sig, group=Sample), position = position_dodge2(width = 0.7)) +
  scale_fill_gradient(low = "#2196F3", high = "#FF5722", guide = "none") +
  scale_color_gradient(low = "#2196F3", high = "#FF5722", guide = "none") +
  scale_alpha_discrete(range = c(0.2, 0.7), guide = "none") +
  scale_y_continuous(expand = c(0, 0), breaks = c(-0.3, 0, 0.3)) +
  facet_grid(Component ~ var_Interoception, scales = "free_y", switch = "y") +
  coord_flip(ylim = c(-0.45, 0.45)) +
  theme_modern(axis.title.space = 5) +
  theme(panel.grid.major.y = element_line(),
        axis.title.y = element_blank(),
        axis.line.x = element_blank(),
        axis.text.y = element_text(size = rel(0.5)),
        axis.text.x = element_text(size = rel(0.8)),
        plot.title = element_text(face = "bold", hjust = 0),
        strip.placement = "outside",
        strip.background.y = element_rect(fill = "#E0E0E0", color = "white"),
        strip.background.x = element_rect(fill = "#F8BBD0", color = "white")) 

Warning: Using alpha for a discrete variable is not advised.

ggsave("figures/Figure2.png", p, width=21*0.6, height=29.7*0.45, dpi=200, bg="white")
p

Discussion

Benefits of the IAS: - Straightforward and sensation-centered items

Recommendations: - Remove Itch (redundant + issue in Chinese) - Use analog scales

Limitations: - Not much clear theorethical or empirical structure (small grouping of items) - Limited variability (clear mode at 4/5) - Ambiguous items which grouping depends on the context (and its awareness) - E.g., heart beating fast, vomit when nauseaous - Few items for some modalities (e.g., 1 for heart) - Positive phrasing: benefits but also might exacerbate positivity bias (and thus unidimensionality)

Need for context-specific items (cross-modal when possible, i.E., cardioception, respiroception, etc.).

New Scale: Multimodal Interoceptive Sensitivity Scale (MISS)